function uy_d = uy_d(x,out_mu)
%%%%%%%%%%true solution%%%%%%%%
epsilon=1E-7;
 n=size(x,1);
 r=0.5;
uy_d=zeros(n,1);
 for i=1:n
     if  x(i,1)^2+x(i,2)^2 - r^2>epsilon
       uy_d(i,1)=2*(x(i,1)^2+x(i,2)^2)/out_mu+(1-2/out_mu)*r^2+ ...
                log(1+x(i,1)^2+x(i,2)^2)/out_mu+(1-1/out_mu)*log(1+r^2);
     else 
         uy_d(i,1) = x(i,1)^2+x(i,2)^2+log(1+ x(i,1)^2+x(i,2)^2);
     end
         

end

         



end

